library(ruta)
library(keras)
library(rARPACK)
library(ggplot2)
#library(tensorflow)
#install_tensorflow()
### synthetic data
load("data_AE_NL.RData")
colnames(data.yoni)

#x1 <- seq(0,1,length=100)
#x2 <- x1**2 -1

x1 <- seq(-2,2,length=100)
x2 <- tanh(x1)
x1 <- x1/2
x1 <- x1+ rnorm(100,sd=0.02)
x2 <- x2+rnorm(100,sd=0.02)
x1 <- x1/max(range(x1)) 
x2 <- x2/max(range(x2)) 

#X1 <- (x1+x2)/sqrt(2)

#X2 <- (x2-x1)/sqrt(2)
### After rotation 
data.yoni <-read.csv("data_NEW_AE_NL.csv")
x2 <- (data.yoni[,1]+data.yoni[,2])/sqrt(2)
x1 <- (data.yoni[,2]-data.yoni[,1])/sqrt(2)
#data.gen <- read.csv(file = "2d_data_forAE.csv")
x2 <- x2/max(range(x2)) 
#X1 <- x2
#X2 <- x1
data.gen <- data.frame(x=x1,y=x2)
x1<- c(-0.06546971878458913, 0.08955351424364644, 0.16665722052141824, 0.11411657942376874, 0.1363071321728137, 0.1955090740698739, 0.31084996162819656, 0.14231549753234274, 0.26544270557038946, 0.47022857754133096, 0.41715561073563884, 0.4142167125902469, 0.4520941716362304, 0.44932857799543163, 0.47241161316824704, 0.5321876563889443, 0.44131206626845354, 0.6714078541000028, 0.4906215936733046, 0.6707004248895789, 0.5270382003022934, 0.6928687276823339, 0.7101584950352954, 0.6985078767079299, 0.6383100048261006, 0.6967320747364463, 0.6744612090900675, 0.733911082036169, 0.6568714963148483, 0.5965717997118055, 0.7650393809330761, 0.7114871521542336, 0.6535809748628225, 0.6208124172266977, 0.6680478099736776, 0.5364643863836802, 0.5927352797383787, 0.5367040789303685, 0.5344018853568318, 0.4766215687099573, 0.34890231408310673, 0.4527912503705071, 0.3914588451256399, 0.19631448817325153, 0.315410895087687, 0.20069888861626295, 0.12199803537478851, 0.1412079645477069, 0.11368539729386576, 0.022077085489203833, -0.12986968425140777, -0.13607278768348635, -0.11033814428209501, -0.1615077471100282, -0.21592640049645623, -0.21389664449547052, -0.3403552374125211, -0.3571571808744574, -0.3156502497952868, -0.49076983995953194, -0.4820587816494964, -0.6197576644589067, -0.619047014573977, -0.5564304886728231, -0.4402241435152499, -0.6543618379392079, -0.5316412913959869, -0.7508002947581249, -0.6984312902591951, -0.6511547031338045, -0.7576162536306754, -0.6878001497568512, -0.7216653531599665, -0.7269622090088661, -0.6998473641066769, -0.6252269165561817, -0.6306359032344657, -0.6987656587615334, -0.6030030369274454, -0.608207246635394, -0.5968267021599665, -0.6018578811788621, -0.6668424586493785, -0.5896786067539002, -0.5214313980300975, -0.467748896462143, -0.5718899508703494, -0.47747679252446834, -0.47952284918438715, -0.5904146838250022, -0.33545604162451337, -0.2575910336962513, -0.29090286412431426, -0.23873383932483377, -0.30125049124795705, -0.18365212211762744, -0.1859789143426045, -0.07729802441542394, 0.08276999063370333, -0.016559214840358105)
x2 <- c(-0.8171667942259502, -0.8198440970049942, -0.8175392915717841, -0.8104204814567585, -0.7759709015927896, -0.7890928732234493, -0.7545468020821812, -0.7835298544997469, -0.7535329323145279, -0.7374356655254856, -0.7306736267324526, -0.7285882300492885, -0.7000222784289601, -0.6783092169817261, -0.6742501646120124, -0.6562701204027782, -0.6501432144761741, -0.651334928169679, -0.6303591209670334, -0.6273260190899269, -0.6128613497996575, -0.5859053168653625, -0.5776279054015778, -0.5744601193751447, -0.5624634352842124, -0.5404957277496986, -0.5020357762406819, -0.5227555148694714, -0.47371198355270033, -0.4681779911310239, -0.44343658523260493, -0.43939972376378417, -0.4290257965096936, -0.4036687958165743, -0.3905634781054041, -0.35010307817166575, -0.3330469252674355, -0.31824622103523015, -0.28146090460239986, -0.2802271365457384, -0.23268340348690345, -0.19421128645071015, -0.20701636163047546, -0.16497132525265315, -0.15236837210348, -0.10759748874779197, -0.07607467111210785, -0.0536742094456272, -0.04026252595442332, -0.009673686644587223, 0.04474325783333033, 0.04680422574353147, 0.0664029581065565, 0.11814315300457073, 0.14685683445655207, 0.11804928477588915, 0.15622339374807767, 0.17709545937755217, 0.2274552054952611, 0.2266936321978855, 0.2701810669866687, 0.29243418910720687, 0.2953428833596076, 0.33397012683304445, 0.36756553441457496, 0.3724258034937967, 0.39033615449545656, 0.4353356339642881, 0.4128055320094032, 0.44418911070722134, 0.45440726376081214, 0.47207141766645966, 0.4883403592208842, 0.531498482306488, 0.5176381293767588, 0.5577473010384588, 0.5456940992922853, 0.5497243364753147, 0.6001255527612617, 0.6020086923498202, 0.6068088531293669, 0.6325013025836388, 0.6383668537347054, 0.6550963362457112, 0.6651074856365792, 0.6705890275446231, 0.7070497543004073, 0.701368552955875, 0.7115930644355478, 0.7181903985977272, 0.7257975662517651, 0.7382078224435049, 0.7495045101949162, 0.7720212590602913, 0.7968851960670642, 0.7876304691009789, 0.8161879515357576, 0.8200049539641215, 0.834338996172365, 0.828795980483226)
### PCA
mus = colMeans(data.gen)
x_train_c =  as.matrix(sweep(data.gen, 2, mus))
SVDS = svd(x_train_c)

ZpcaTrain = x_train_c %*% SVDS$v[,1] # PCA projection of test data
Rpca = sweep(ZpcaTrain %*% t(SVDS$v[,1]), 2, -mus) # PCA reconstruction
dim(Rpca)

error.reconstr.PCA <- sum((Rpca-data.gen)^2)
#data.gen1 <- data.frame(data.gen,color=rep("black",200))
#colnames(data.gen1)[1:2] <- c("x","y")
dat.PCA <- data.frame(PC1=Rpca[,1],PC2=Rpca[,2])


### AE
network <- input() + dense(1, "tanh") + output("tanh")

network.simple <- autoencoder(network)#, loss = "binary_crossentropy")
model = ruta::train(network.simple, as.matrix(data.gen,ncol=2,nrow=dim(data.gen)[1]), epochs = 2500)
decoded.simple <- reconstruct(model,as.matrix(data.gen,ncol=2,nrow=dim(data.gen)[1]))
colnames(decoded.simple) <- c("x1","x2")
decoded.simple <- as.data.frame(decoded.simple)
error.reconstr.AES <-sum((decoded.simple-data.gen)**2)

## get the code
res.code.AE <- encode(model,as.matrix(data.gen,ncol=2,nrow=dim(data.gen)[1]))


error.reconstr.PCA <- sum((Rpca-data.gen)^2)
dat.PCA <- data.frame(PC1=Rpca[,1],PC2=Rpca[,2])


p <- ggplot()+
  geom_point(data = dat.PCA, aes(x = PC1, y = PC2, colour="PCA"),size = 1.5)+
  #scale_fill_manual(name = "", values = c("PCA" = "red"))+
  geom_point(alpha = 0.8,) +
  geom_point(data=data.gen,aes(x = x,y = y,fill="Data points"),colour="black", size = 1.5)+
  # scale_color_manual(values = c("black" = "black"))+
  #geom_point(data=data.true,aes(x = x,y = y,colour="True pattern"), size = 1.5)+
  geom_point(data=decoded.simple,aes(x = x2,y = x1,colour="Shallow AE"), size = 1.5)+
 # geom_point(data=decoded2,aes(x = x1,y = x2,colour="Deep AE"), size = 1)+
  labs(x = expression(x[1]), y=expression(x[2]))+
  theme_bw() +
  theme(axis.text = element_text(size = 16),axis.title  = element_text(size = 16),
        legend.title =element_blank(),
        legend.text = element_text( color = "Black", size = 16),
        legend.position = c(0.15,0.85))

p


ggsave("../../../figures/chapter_3_figures/manifold_synthetic_3.pdf", width = 8, height = 8)
ggsave("manifold_synthetic_rotate.pdf", width = 8, height = 8)


data.yoni <- cbind(data.gen[,1],data.gen[,2],Rpca[,1],PC2=Rpca[,2],decoded.simple$x1,decoded.simple$x2,res.code.AE,ZpcaTrain)
colnames(data.yoni) <- c("x1","x2","x1hat.pca","x2hat.pca","x1hat.AE","x2hat.AE","code.AE","code.PCA")
save(data.yoni,file="data_AE_NL.RData")
write.csv(data.yoni,file="data_NEW_AE_NL.csv",row.names = FALSE)
